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ABSTRACT 

We investigate the migration of massive extrasolar planets due to gravitational inter- 
action with a viscous protoplanetary disc. We show that a model in which planets 
form at 5 AU at a constant rate, before migrating, leads to a predicted distribution 
of planets that is a steeply rising function of log(o), where a is the orbital radius. 
Between 1 AU and 3 AU, the expected number of planets per logarithmic interval in 
a roughly doubles. We demonstrate that, once selection effects are accounted for, this 
is consistent with current data, and then extrapolate the observed planet fraction to 
masses and radii that are inaccessible to current observations. In total, about 15% of 
stars targeted by existing radial velocity searches are predicted to possess planets with 
masses 0.3Mj < M p sm(i) < lOMj , and radii 0.1AU < a < 5AU. A third of these 
planets (around 5% of the target stars) lie at the radii most amenable to detection via 
microlensing. A further 5-10% of stars could have planets at radii of 5AU < a < 8AU 
that have migrated outwards. We discuss the probability of forming a system (akin to 
the Solar System) in which significant radial migration of the most massive planet does 
not occur. About 10-15% of systems with a surviving massive planet are estimated to 
fall into this class. Finally, we note that a smaller fraction of low mass planets than 
high mass planets is expected to survive without being consumed by the star. The 
initial mass function for planets is thus predicted to rise more steeply towards small 
masses than the observed mass function. 

Key words: accretion, accretion discs — solar system: formation — planetary sys- 
tems: formation — planetary systems: protoplanetary discs — gravitational lensing 



1 INTRODUCTION 

Radial velocity surveys of nearby stars show that a signif- 
icant fraction - at least 8% - have massive planets with 
orbital radii substantially less than that of Jupiter (Marcy 
& Butler 2000; Udry et al. 2000; Butler et al. 2001). Most 
of these planets lie at radii where massive planet forma- 
tion is theoretically believed to be difficult (Bodenheimer, 
Hubickyj & Lissauer 2000). The difficulties are most pro- 
nounced for those planets at the smallest orbital radii, where 
the temperatures in the protoplanetary disc would have ex- 
ceeded those for which ices (and possibly even dust) can 
exist. These problems can be avoided if planets formed at 
larger radii, and then migrated inwards to where they are 
detected today. Mechanisms that could lead to this migra- 
tion include gravitational interaction with a gaseous and vis- 
cous protoplanetary disc (Lin, Bodenheimer & Richardson 
1996), planet-planet scattering (Rasio & Ford 1996; Weiden- 
schilling & Marzari 1996; Lin & Ida 1997; Ford, Havlickova 



& Rasio 2001; Papaloizou & Terquem 2001), or scattering 
of planetesimals by massive planets (Murray et al. 1998). 

In the right circumstances, all of the suggested mech- 
anisms can lead to substantial migration. More quantita- 
tive comparisons with the data are therefore warranted. Re- 
cently, Trilling, Lunine & Benz (2002) have reported promis- 
ing results for the protoplanetary disc migration model. 
They showed that the broad distribution of orbital radii 
was consistent with the planets becoming stranded, during 
their inward migration, by the dispersal of the protoplane- 
tary disc. They used the model to estimate both the fraction 
of stars that must have formed planets, and how many ought 
to have survived migration. 

In this paper, we extend the work of Trilling et al. 
(2002). We include a physical mechanism for disc disper- 
sal into a model for the formation and migration of mas- 
sive planets, and make quantitative comparison with the 
observed distribution of planetary orbital radii. Our disc 
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model, described in Section 2, combines viscous evolution 
with mass loss at the outer edge, for example due to photo- 
evaporation (Johnstone, Hollenbach & Bally 1998; Clarke, 
Gendrin & Sotomayor 2001). We then run this model re- 
peatedly, on each occasion adding a single planet to the disc 
at a specified formation time ti OTm . By varying ti otm , sub- 
ject to the constraint that no planets can be formed once 
the disc mass is too low, we study in Section 3 how the fi- 
nal orbital radii of planets depend upon the formation time. 
With the addition of plausible assumptions about how the 
rate of planet formation varies with time, we then convert 
the results into a prediction for the radial distribution of 
planetary orbits. In Section 4 we compare this distribution 
with the current data, and find that good agreement is ob- 
tained. In Section 5 we make use of the model to estimate 
the fraction of stars targeted by radial velocity surveys that 
ought to harbour planets, including those with masses too 
low, or orbital radii too large, to be currently detected. We 
also predict the number of stars with planets at radii suitable 
for detection via microlensing, and the fraction of systems 
where significant planetary migration does not occur. Our 
conclusions are summarized in Section 6. 



2 METHODS 

2.1 Fully viscous disc model 

The simplest protoplanetary disc model that meets the basic 
observational constraints on disc mass, disc lifetime, and disc 
destruction time-scale is a viscous disc with mass loss from 
the outer regions (Hollenbach, Yorke & Johnstone 2000; 
Clarke, Gendrin & Sotomayor 2001). We use this model for 
the majority of the calculations described in this paper, since 
it has a minimal number of free parameters. 

For the viscous disc model runs, we assume that over 
the relevant range of radii (roughly, between 0.1 AU and 
10 AU), the kinematic viscosity v can be described as a 
power law in radius R, 



v = Vq 



(l Au) 



(1) 



The initial conditions for the surface density are a constant 
accretion rate disc, with a profile which corresponds to the 
assumed viscosity of equation (|lj). At radii R ^> i£j n) the 
inner edge of the disc, this implies, 



M y R-e 
- — = L.oR 

i-KV 



(2) 



where M is the accretion rate and Eo a constant. For our 
standard disc model, we take j3 = 3/2 (Weidenschilling 1977; 
Hayashi 1981), and choose vq = 1.75 x 10 13 cm 2 s -1 to give 
a sensible time-scale for disc evolution of a few Myr. Fi- 
nally, Eo is set in the initial conditions so that the initial 
disc mass is O.IMq. The inner and outer boundaries for the 
calculation are Ri n = 0.066 AU, and R ou t = 33.3 AU. To 
extend the predictions to smaller radii, a limited number 
of runs have also been completed with an inner boundary 
at Ri n — 0.03 AU. A zero torque boundary condition (i.e. 
a surface density E = 0) is applied at the inner boundary, 
while at the outer boundary we set the radial velocity to 
zero. 



Current observations do not directly constrain the sur- 
face density profile of protoplanetary discs at the AU scale. 
To check how sensitive our results are to changes in this 
profile, we have also run a model with /3 — 1/2, again nor- 
malized to give a sensible time-scale for disc evolution. A 
flatter profile than the R~ 3 ^ 2 of our standard model is ex- 
pected if the angular momentum transport within the disc 
follows the Shakura & Sunyaev (1973) a-prescription, with 
a constant a (e.g. Bell et al. 1997). 

Even a small mass of gas, trapped exterior to the orbit 
of a massive planet, will eventually soak up a large amount of 
angular momentum and drive migration. This gas has to be 
lost in order to obtain a converged final radius for migrating 
planets. Trilling et al. (2002) make the simplest assumption, 
and halt migration by removing the disc instantaneously. We 
adopt a different approach, and assume that the gas is re- 
moved by photoevaporation. In this process, the surfaces of 
the disc are heated by the absorption of ultraviolet radiation, 
which may originate either from the central star or from an 
external source (e.g. a massive star in a young cluster). The 
hot gas can escape as a wind at radii R > R g , where the 
escape velocity is less than the sound speed of the heated 
material. Observations of the Orion nebula show that this 
process drives mass loss at significant rates, at least in clus- 
ters containing massive stars (Bally, O'Dell & McCaughrean 
2000). 

Our implementation of photoevaporation within a time- 
dependent disc model is based upon the description given 
by Johnstone, Hollenbach & Bally (1998; see also Shu, John- 
stone & Hollenbach 1993; Richling & Yorke 2000). The stan- 
dard model has R g = 5 AU, and a radial distribution of mass 
loss, 



RT 



R > Rn 



E = R< R„ 



(3) 



For a disc with an outer radius U > fl 9 , this means that 
the total mass loss scales linearly with R. The initial mass 
loss rate (integrated over the disc out to 25 AU) is chosen 
to be M = 5 X 10 _9 M Q yr _1 . 

Figure [j] shows how the disc mass changes with time, 
for the standard model with (3 = 3/2 and mass loss from 
R > 5AU. The disc mass declines from 0.1M Q to 10" 2 M o in 
3 Myr, mainly due to accretion, before dropping rapidly once 
the effects of the disc wind take hold. Although the power 
laws and location of the turnover vary with the adopted 
parameters, generically similar behavior occurs for a range 
of disc viscosity and mass loss rates (Clarke, Gendrin & 
Sotomayor 2001). 



2.2 Magnetically layered disc model 

In order to check the robustness of the results, we investigate 
how the results would change if we adopted a fundamen- 
tally different, but strongly theoretically motivated, model 
for disc evolution. The magnetically layered model proposed 
by Gammie (1996) is based upon the observation that cold 
gas at T £ 10 3 K has a low ionization fraction, which sup- 
presses angular momentum transport via magnetohydrody- 
namic turbulence (Fleming, Stone & Hawley 2000). If there 
are no other sources of angular momentum transport, then 
only the hot inner disc, plus a surface layer ionized by cosmic 
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Figure 1. Evolution of the mass of the standard disc model with 
time. Up to around 3 Myr, the disc mass decays as gas is accreted 
onto the star, while the effects of the disc wind remain small. At 
later times, a steeper decline occurs as mass loss truncates the 
outer radii of the disc. The shaded band shows the range of time 
over which a 2Mj planet can be formed at 5 AU and survive with 
a final orbital semi-major axis 0.1AU < a < 6AU. 

rays, will be viscous. At radii R ~ 1AU, this implies that ac- 
cretion occurs only through a thin surface layer, below which 
lies a thick 'dead zone' of quiescent gas. Numerical models 
for the long-term disc evolution (Armitage, Livio & Pringle 
2001) find that there is an early phase, lasting for perhaps a 
Myr, during which accretion occurs in short high accretion 
rate outbursts. This is followed by a quiescent phase, char- 
acterised by low rates of accretion onto the star. This is due 
to a much reduced (compared to fully viscous disc models) 
viscosity at radii of the order of 1 AU. 

Detailed models for the layered disc are able to pro- 
vide a reasonable fit to observations of the accretion rate in 
Classical T Tauri stars (Gammie 1996; but see also Stepinski 
1998). For the purpose of studying planetary migration, we 
adopt an approximate approach, and assume that the most 
important effect of the layered disc is to reduce the viscos- 
ity near the midplane. Specifically, after 1 Myr, we reduce 
the viscosity where the surface density exceeds the thickness 
of the magnetically active surface layer. The new viscosity 

^layer IS, 

. / lOOgcm" 2 A fA . 

flayer = V X mm I , 1 I (4) 

where v is as given in equation (Q), and the same viscosity 
parameters as for the standard disc model are chosen. To 
obtain a sensible disc lifetime we also increase the mass loss 
rate to M — 3 x 10 -8 M©yr -1 , and allow mass to be lost at 
all radii by taking R g = R[ n . 



2.3 Planet formation assumptions 

The location in the disc where massive planet formation is 
easiest depends upon two competing factors. The character- 
istic time-scales for planet building in the disc are faster at 
small radii. However, most of the mass in protoplanetary 
discs lies at large radii, with a jump in the surface density 
of solid material beyond the 'snow line' where ices can form 
(Hayashi 1981; Sasselov & Lecar 2000). We assume that the 
outcome of this competition is that planets form at an or- 
bital radius a — 5AU, similar to that of Jupiter in the Solar 
System. The exact choice of this location is less important 
than the fact that it lies outside the radii where extrasolar 
planets are currently observed. In the model which we are 
testing, this means that migration is necessary to explain 
the orbital radii of all currently known planets. 

In each run of the disc model, we form a planet with 
mass M p , at time t[ OTrn (where t — is defined as the time 
when the disc mass is O.IM©). We specify the mass of our 
planets (typically 2Mj, where Mj is the mass of Jupiter) in 
advance of each run, but test at the time of formation that 
there is enough mass in the vicinity of the planet to form 
it consistently. Specifically, we check that the disc mass be- 
tween R = 0.6a and R = 1.6a exceeds the desired planet 
mass. Provided that this is satisfied, we remove the appro- 
priate amount of gas from the disc around R = a, and 
add a planet with mass M p at that location. This instan- 
taneous formation scheme ignores the possibly lengthy time 
required to assemble the planet's gaseous envelope (e.g. Lis- 
sauer 1993, and references therein). However, Trilling, Lu- 
nine & Benz (2002) have shown that including a gradual 
build-up of mass makes little difference to the outcome of 
migration. 

After running the model a number of times, the basic 
output is a plot of the final radius of a planet as a function 
of the time when it formed. To convert this to the number of 
planets expected at different radii, which is the observable 
quantity, we also need to specify how the probability of mas- 
sive planet formation, per unit time, varies with time. We 
assume that this probability is uniform (or equivalently, that 
the rate of massive planet formation, averaged over many 
stars, is constant). This cannot be true over long periods 
of time. However, we are only interested in the window of 
formation times which allow a massive planet to survive mi- 

S-ation. This window, shown as the shaded region in Figure 
is short compared to the disc lifetime, and as a result there 
is only a relatively small change in the physical properties of 
the disc during this time. For example the disc mass, which 
drops by two orders of magnitude over the first 3.5 Myr, 
changes by less than 50%. Whatever factors influence the 
probability of planet formation, it should therefore be a rea- 
sonable first approximation to assume that the probability 
is constant over the short range of formation times which 
lead to a surviving massive planet. 

2.4 Model for planetary migration 

Our model for planetary migration is almost identical to that 
of Trilling et al. (2002), who describe in detail the approx- 
imations involved. The coupled evolution of the planet and 
the protoplanetary disc is described by the equation (Lin & 
Papaloizou 1986), 
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dE _ 1 d 
~dt ~ RdR 



2ASfl 3/2 

(GM„y/ 2 



(•») 



The first term on the right-hand side of this equation de- 
scribes the diffusive evolution of the surface density due to 
internal viscous torques (e.g. Pringle 1981). The second term 
describes how the disc responds to the planetary torque, 
A(R,a), where this function is the rate of angular momen- 
tum transfer per unit mass from the planet to the disc. We 
take, 



A = 



A = 



q 2 GM„ ( R 



2R 

q 2 GM, 
2R 



A" R < a 



aJ R>a 



where q = M p /M*, 
A p = max(_ff, \R — a\), 



(6) 



(7) 



and H is the scale height of the disc. Guided by detailed 
protoplanetary disc models (Bell et al. 1997), we adopt H = 
0.05ii. This form for A (equation ^ is that used by Lin & 
Papaloizou (1986), modified to give a symmetric treatment 
for the disc inside and outside the orbit of the planet. 

The transfer of angular momentum leads to orbital mi- 
gration of the planet at a rate, 



do 
dt 



GAL 



1/2 / 47r 



- | / *AEd*. 



(8) 



Although simple, this formalism for treating the planet-disc 
interaction has been shown (Trilling et al. 1998) to give com- 
parable results to more sophisticated methods (Takeuchi, 
Miyama & Lin 1996). 

We solve equation Q5j) using an explicit method on a grid 
that is uniform in a scaled variable oc _R 1//2 , with 300 grid 
points between R ln and -Rout - The stellar mass is M» = Mq. 
Typically, the timestep is limited by the radial velocity in 
the disc gas very near the location of the planet. Evolving 
the system on this (small) timestep serves no useful purpose 
once a gap has been opened, so we reduce A, and thereby 
limit |u r |, in the vicinity of the planet. The results of the 
migration calculations are unaffected by this modification. 



3 RESULTS 

3.1 Predicted semi-major axis of stranded planets 

Figure |^ shows the final semi-major axis of planets as a func- 
tion of the formation time, using the standard disc model of 
Fig. [lj. Each point represents a separate run of the model. 
The fate of a planet depends upon when it was formed, and 
upon its mass. Planets that are formed too early migrate 
inward of 0.1 AU, the smallest radius we consider quanti- 
tatively in this study. What happens to them subsequently 
depends upon whether there is a stopping mechanism to halt 
further migration, but many of them seem likely to be con- 
sumed by the star. Planets that start to form too late can 
migrate outwards, because at late times the gas in the disc at 
5 AU is itself moving away from the star. At still later times, 
it becomes impossible to form a Jupiter mass planet at all, 
because the dispersing disc has too little surface density (c.f. 
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Figure 2. The final semi-major axis of planets, following migra- 
tion, as a function of the planet formation time (i.e. the epoch 
when the planet formed at 5 AU). The standard disc model with 
(3 = 3/2 and mass loss at R > 5AU was used. The solid curve 
shows results for 2Mj planets, the dashed curve results for lMj 
planets. More massive planets can form earlier without being con- 
sumed by the star. 



the discussion in Shu, Johnstone & Hollenbach 1993). In be- 
tween these limits, there is a window of formation times, 
lasting for around 0.5 Myr for a 2Mj planet, which result in 
the planet being stranded at radii between 0.1 AU and 6 AU. 
The migration of more massive planets slows down earlier, 
since their larger angular momentum produces greater re- 
sistance to disc- induced migration (Syer & Clarke 1995). 
Hence, the window of allowable formation times is longer 
for more massive planets than for less massive ones. 

3.2 Dependence upon disc model and planet mass 

Figure ^ shows how the final planetary radius depends upon 
the formation time (note that we plot a(if orm ) on a logarith- 
mic scale in this figure). The curve steepens towards small 
radii / early formation times. This implies that if planet for- 
mation occurs with uniform probability per unit time during 
the allowed window, then fewer planets per logarithmic in- 
terval in a are expected at small radius (i.e. equal intervals 
in log(a) correspond to smaller intervals in Aif orm at small 
a). 

With the assumption of a uniform rate of planet forma- 
tion, the predicted number of planets per logarithmic inter- 
val in radius is, 



diV D 



dlog(a) 



(9) 



dlog(a) y dt form 

which can readily be derived from Fig. |j| Figure |^ shows the 
results for two different planet masses, and for two variants 
on the standard disc model. A cubic spline fit has been used 
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Figure 3. Predicted number of extrasolar planets per logarith- 
mic interval in semi-major axis. These curves are derived from the 
a(tf orm ) results shown in Fig. 3 by making the additional assump- 
tion that the rate of planet formation is constant. The solid curve 
shows results for the standard disc model (2Mj planets, f) = 3/2, 
R g = 5AU). The remaining curves show the effect of changes to 
these parameters: a lower planet mass M p = lMj (short dashes), 
increased radius of mass loss R g = 10AU (long dashes), and dif- 
ferent viscosity f3 = 1/2 (dot-dashed). Because the (3 = 1/2 model 
has a lower surface density at 5 AU, the planet mass used was 
reduced to 0.5Mj. The absolute number of planets is arbitrary, 
and has been normalised to unity at 1 AU. 



to obtain a smooth estimate of the derivative in equation 
(^), and the curves have been normalised to unity at 1 AU. 

The predicted number of planets per logarithmic inter- 
val in a rises rapidly with increasing radius for all the models 
considered. For the standard disc model, the predicted num- 
ber rises by more than a factor of two between 0.1 AU and 
1 AU, and by a further factor of two by 3 AU. The good 
agreement between the curves for lMj and 2Mj planets 
shows that the observed mass function for massive planets is 
predicted to be the same at different radii. This means that 
we can use the observed mass function at small radii, which 
is complete down to lower masses, to estimate how many low 
mass planets must accompany known high mass planets at 
larger radii. It does not mean that the initial mass function 
for planets is the same as the observed one. In fact, Fig. ^ 
makes it clear that high mass planets can form over a longer 
interval than low mass planets and still survive migration. 
To end up with the same number of lMj as 2Mj planets 
per unit interval in log(a), more low mass planets must have 
formed. In other words, the (currently unobservable) initial 
mass function for planets must rise more steeply towards 
low masses than the observed one. 

Figure ^ also shows the predicted relative abundance of 
planets for two variants of the disc model. Increasing the in- 
ner radius beyond which mass is lost to 10 AU does not alter 




a / AU 



Figure 4. Predicted number of extrasolar planets per logarith- 
mic interval in semi-major axis. The solid line shows the results 
for the standard disc model, as in Fig. 3. The short-dashed and 
long-dashed curves show models in which the disc is lost instan- 
taneously at t = 3 Myr, or t = 4 Myr, respectively. The predicted 
distribution of planets that have migrated inwards is almost the 
same in all of these models. 

the predicted distribution of planets within 5 AU. Likewise, 
the predicted distribution is not substantially changed in a 
model where the power-law slope of the viscosity is /3 = 1/2. 
Note that for the /3 — 1/2 run we considered planets with 
mass 0.5 Mj, since insufficient gas was available between 
0.6a and 1.6a to form planets of several Jupiter masses. 



3.3 Comparison with an instantaneous disc 
dispersal model 

In Figure ||, we show how the model including disc mass loss 
compares to a simpler model in which the disc is assumed 
to be lost instantaneously at a time tdisperse- The disc model 
used is identical to our standard model with j3 = 3/2, except 
that the rate of mass loss is set to zero. We compute two 
models, with tdispcrsc equal to 3 Myr and 4 Myr respectively. 

From Fig. kl it can be seen that the instantaneous disc 
dispersal model makes very similar predictions for the ra- 
dial distribution of inwardly migrating planets. Larger differ- 
ences, however, occur for outward migration. In the instan- 
taneous disc dispersal model, there is no outward migration 
of planets, because the radial velocity in the disc at 5 AU 
remains inward at tdispcrsc- The model including mass loss, 
conversely, allows a substantial fraction of planets to migrate 
outwards. From Fig. ^, we find that for Jupiter mass planets 
as many as a third may end up in orbits exterior to that in 
which they formed. 
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Figure 5. Results of planetary migration in a magnetically lay- 
ered disc, where we have assumed that the layer is established at 
t = 1 Myr. The upper panel shows the disc mass as a function of 
time, with the shaded band illustrating the range of planet for- 
mation times that result in the planet having a final semi-major 
axis between 0.1 AU and 6 AU. The lower panel shows the final 
semi-major axis as a function of formation time. 



3.4 Effect of ongoing planetary accretion 

Even after a gap has been opened, numerical simulations 
show that planetary accretion may continue (Bryden at al. 
1999; Kley 1999; Lubow, Seibert, & Artymowicz 1999; Nel- 
son et al. 2000). Significant ongoing accretion - which is not 
included in our models - would lead to a correlation between 
planet mass and orbital radius, with more massive planets 
expected at small orbital radius. There is no evidence for 
such a correlation in the data, which could either indicate 
that migration has not occurred, or that ongoing accretion 
occurred only at a low level. 

Lubow, Seibert, & Artymowicz (1999) find that, for low 
masses, the rate of ongoing accretion can be comparable to 
the disc accretion rate outside the planet. This is of the 
order of 10 -9 M0yr _1 at the epoch when surviving planets 
form. Migration from 5 AU takes a few hundred thousand 
years, so we estimate that accretion across the gap could 
add a few tenths of a Jupiter mass to low mass planets. 
This is a relatively small fraction of the mass for the 2Mj 
planets that are the main focus of this paper. However, it is 
only marginally consistent with the existence of the lowest 
mass hot Jupiters, which have M v sin(i) < 0.3 Mj. This 
may indicate that these planets formed within the snow line 
(leading to a shorter migration time). 



3.5 Magnetically layered protoplanetary discs 

Figure ^ shows the results for the runs using the magneti- 
cally layered disc model. The decline of the disc mass with 



time shows the effect of the switch to a low viscosity state. 
After the initial fully viscous phase, which resembles the 
standard model, there is a long plateau which arises because 
the low viscosity in the inner disc produces a bottleneck to 
accretion onto the star. During this period the disc is slowly 
destroyed from the outside inwards by the disc wind, which 
finally disperses the disc completely after 6 Myr. We note 
that, in a magnetically layered disc model, high rates of 
mass loss from the disc are needed in order to disperse the 
disc within a reasonable time-scale. The mass loss rates in- 
ferred for discs in parts of the Orion nebula are easily large 
enough (Johnstone, Hollenbach & Bally 1998), but it is un- 
known whether high enough rates are possible for isolated 
discs. 

Figure ^ also shows the final location of migrating plan- 
ets. As before, we assume these planets all formed at 5 AU, 
and vary the time of formation in a series of runs. The result- 
ing a(tf orm ) curve flattens for a > 1AU. As with the standard 
disc model, therefore, we predict a broad distribution of fi- 
nal orbital radii, with more planets (per logarithmic interval 
in a) ending up at large radii than at small radii. However, 
there are also important differences. Planets must form ear- 
lier (relative to the disc lifetime) in the layered disc if they 
are to migrate successfully to small orbital radii. This is 
because little migration can take place once the disc has en- 
tered its quiescent, low viscosity phase. For a planet to end 
up inside 1 AU, it must have formed at an earlier time when 
the disc was still viscous and able to drive rapid migration. 



4 COMPARISON WITH CURRENT DATA 

Before comparing our results to current data, we need to 
consider possible biases. In an idealized radial velocity search 
for extrasolar planets, two important selection effects are ex- 
pected. First, there is almost no sensitivity to planets with 
orbital periods P longer than the duration of the survey (un- 
til a whole orbit is seen, the planetary signal may be con- 
fused with that from a more massive, distant companion). 
Second, a planet with mass M p and orbit inclination i pro- 
duces a radial velocity signal with an amplitude proportional 



to M v sin(i)a 



-1/2 



A survey which makes a fixed number of 



observations, with some given sensitivity, can therefore de- 
tect planets above a minimum M v sin(i) cx a 1 ' 2 . Although 
this is a relatively weak radial dependence, it still corre- 
sponds to an order of magnitude difference in the minimum 
detectable mass over the range of radii probed by current 
surveys. 

Real surveys are less easily characterised. Instruments 
have improved steadily over time, and the intrinsic limits to 
radial velocity measurements vary on a star- by-star basis. 
Nevertheless, the selection effect caused by the scaling of 
the minimum detectable mass with a 1//2 can readily be seen 
in the data. Figure ^ shows the location of all the known 
extrasolar planets with M v sin(i) < 10 Mj, in the M p sin(i) 
- a plane. Most of the hot Jupiters at a < 0.1 AU, where 
surveys are most sensitive, have masses that are actually 
substantially smaller than that of Jupiter, while the typical 
detected planet with an orbital radius greater than 1 AU 
has a mass of several Mj . 

A simple count of the total number of detected planets 
as a function of radius leads to the cumulative distribution 
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Figure 6. The observed distribution of extrasolar planets in the 
M p sin(i) - a plane. The diagonal lines indicate the approximate 
selection limits of current radial velocity surveys. The upper line, 
labelled CMB99, shows the minimum mass above which plan- 
ets are detectable at 99% confidence, according to the analysis 
of Cumming, Marcy & Butler (1999). We also plot lines corre- 
sponding to radial velocity amplitudes of 3 ms" 1 and 10.5 ms -1 . 
The rectangular boxes define subsamples of the data, used for 
comparison with the predictions of disc migration models. 



Figure 7. The cumulative probability distribution of extraso- 
lar planets as a function of orbital radius. The dot-dashed curve 
shows the distribution for all detected planets. This distribution 
is consistent (KS test probability P = 0.05) with a uniform dis- 
tribution in log(a), shown as the solid line. We also consider sub- 
samples comprising only those planets with masses above 0.6 Mj 
(long dashes), or 1.2 Mj (short dashes). The hypothesis that these 
subsamples have a distribution that is uniform in log(a) can be 
rejected at more than 99.9% confidence. 



M p sin(i) > 1.2Mj M v sin(i) > 0.6Mj 

0.040 - 3.0 AU 5 x 10 -6 7 x 10~ 4 

0.065- 1.9 AU 10~ 4 10~ 3 



Table 1. Kolmogorov-Smirnov test probability that the observed 
distribution of extrasolar planets is drawn from a distribution uni- 
form in log(a). The subsamples considered have masses exceeding 
0.6 Mj or 1.2 Mj, and orbital radii within the two quoted ranges. 

shown in Figure |^. A Kolmogorov-Smirnov test (e.g. Press et 
al. 1989) shows that this is marginally consistent (KS proba- 
bility P — 0.05) with a uniform distribution in log(a). How- 
ever, from Fig. ^| it can be seen that this uncorrected dis- 
tribution includes low mass planets, at small orbital radius, 
that would not be detectable beyond 1 AU. The mass func- 
tion of planets has substantial numbers of low mass planets 
(Zucker & Mazeh 2001; Tabachnik & Tremaine 2002), so this 
leads to a significant selection effect (Lineweaver & Grether 
2002). 

An unbiased estimate of the radial distribution of ex- 
trasolar planets can be obtained by considering a subsample 
of planets which are sufficiently massive that they could be 
detected at any radius. A conservative definition of such a 
subsample can be obtained by noting that for the Lick planet 
search, circa 1999, the minimum detectable mass at 99% 
confidence is quoted by Cumming, Marcy & Butler (1999) 
as, 



M p sin(i) > 0.7Mj(a/AU) 1/2 . (10) 

Taking this as a guide, a mass cut at M cu t = M p sin(i) > 
1.2 Mj should yield an approximately complete subsample 
out to 3 AU. Since the precision of the observations is likely 
to have improved since 1999, we also consider a subsample 
with the mass cut placed at 0.6 Mj. This corresponds to 
a radial velocity amplitude at 3 AU of 10.5 ms -1 , which 
is equal to the smallest amplitude signal of any detected 
planet. As one might expect, this practical detection thresh- 
old lies at 3 — 4<r, where a ~ 3 ms -1 is the best-case error 
on a single radial velocity measurement. 

As shown in Fig. (?) and Table 1, neither of these subsam- 
ples is consistent with the hypothesis that it is drawn from 
a uniform distribution in log(a), which can be excluded at 
the 99.9% confidence level. This conclusion is also robust to 
different choices of the minimum and maximum orbital radii 
of planets in the subsample. 

To compare with the theoretical distribution derived 
earlier, we henceforth take M cut = 0.6 Mj. This may some- 
what underestimate the number of planets at radii of 2- 
3 AU, but any such bias is smaller than the uncertainties 
introduced by the relatively small sample size. We then 
take the current compilation of 76 detected planets (from 
exoplanets.org), discard those with masses outside the 
range 0.6Mj < M v sin(i) < lOAfj, and bin 54 survivors 
into four logarithmic intervals in semi-major axis, 0.033 AU 
- 0.1 AU, 0.1 AU - 0.3 AU, 0.3 AU - 0.9 AU, and 0.9 AU - 
2.7 AU. 



© 0000 RAS, MNRAS 000, 000-000 



8 P.J. Armitage, M. Livio, S.H. Lubow & J.E. Pringle 




o 1 — 1 — ^ 

0.1 1 
semi-major axis / AU 

Figure 8. Predicted and observed number of extrasolar planets 
as a function of radius. The data (points with error bars showing 
vN errors) show the number of extrasolar planets in four radial 
bins, using a subsample with a mass cut at 0.6 Mj. The theo- 
retical curve is the prediction for planetary migration within the 
standard disc model, with the normalisation chosen to give the 
best fit to the data. 

Figure |§| shows the comparison between observation and 
theory. The disc inside 0.1 AU may well be strongly influ- 
enced by the protostellar magnetosphere (Konigl 1991; Na- 
jita, Edwards, Basri & Carr 2000), which we have chosen 
not to model. We therefore concentrate attention on the 
outer three radial bins, between 0.1 AU and 2.7 AU. The 
number of known extrasolar planets in this (approximately) 
complete subsample rises rapidly with increasing semi-major 
axis. Although the error bars on the data are large, excellent 
agreement is obtained with the predictions of the theoretical 
migration model using the standard disc model. This agree- 
ment persists when different choices of binning and sample 
selection are made. In particular, we can obtain an equally 
good fit (with a different normalization) if we adopt a more 
conservative mass cut at 1.2 Mj. 

From the Figure, we also note that in this mass range, 
the number of hot Jupiters at a < 0.1 AU is in fact consistent 
with an extrapolation of the data from larger radii. There 
is no need for an additional stopping mechanism to explain 
the abundance of these planets at the smallest radii. There 
is, however, weak evidence for an over-abundance of lower 
mass planets (with M p sin(i) < 0.6Mj, not plotted in the 
Figure) in the innermost radial bin. 



5 DISCUSSION 

5.1 Frequency of massive planets 

The good agreement between theory and observations seen 
in Figure @ allows us to use the model to estimate the frac- 



tion of stars targeted by radial velocity surveys that pos- 
sess massive planets, including those with masses too low, 
or semi-major axis too large, to have been detected so far. 
Of course, we cannot extrapolate to masses smaller than 
the ~ Q.2M.J sin(i) that is the smallest yet detected, as we 
have no information at all on the mass function below this 
threshold. We can, however, estimate how many planets with 
0.5Afj sin(i) ought to exist at 1 AU (say), beyond the radius 
where they are currently detectable. 

To do so, we note that the radial velocity surveys are 
most complete (i.e. reach down to the smallest masses) at 
small radii. Between 0.1 AU and 0.2 AU, equation (|l^) sug- 
gests that planets with masses of 0.3Mj sin(i) and above 
are detectable. 8 planets are known in this radial interval, 
corresponding to approximately 1% of target stars (since 
the total haul of 76 planets implies that approximately 
8% of stars have planets). The range of planet formation 
times which result in the planet becoming stranded between 
0.1 AU and 0.2 AU is approximately 2.5 x 10 4 yr, which may 
be compared to the 0.42 Myr window which leaves plan- 
ets between 0.1 AU and 5 AU. Hence, we estimate that 
the fraction of target stars possessing planets with masses 
0.3Mj < M p sin(i) < 10Mj, and radii 0.1AU < a < 5AU, is 
~ 15 %. Roughly half of these have already been detected. 

More speculatively, we can also try to extrapolate to 
larger radii. If we assume that massive planets at 5-8 AU 
formed near the inner boundary of that region and migrated 
outwards, then a similar analysis to that above suggests that 
an additional 5-10% of stars could possess such planets. This 
conclusion is more tentative than for inward migration, how- 
ever, because the extent of outward migration varies signifi- 
cantly with planet mass (Fig. |^) , and could be very different 
in the magnetically layered disc model. 

5.2 Frequency of massive planet formation 

In our standard model, planets which survive migration 
with final orbital radii between 0.1 AU and 6 AU must 
have formed during a 0.5 Myr allowed window. This win- 
dow lasts for only « 20% of the disc lifetime. Migration 
without a stopping mechanism is thus a moderately inef- 
ficient process (Trilling et al. 2002) - several planets are 
likely to have formed for every planet which survives to- 
day. This inefficiency means that either a large fraction of 
stars formed of the order of one planet, or that a smaller 
fraction of stars formed many planets. Our models do not 
distinguish between these possibilities, though in principle 
they make different predictions. For example, the fraction of 
stars which consume planets (and may show resulting metal- 
licity enhancements) is smaller if only a few stars form many 
massive planets. The interpretation of current metallicity 
measurements within such a framework, however, remains 
somewhat controversial (e.g. Gonzalez 1997; Udry, Mayor & 
Queloz 2001; Murray et al. 2001; Pinsonneault, DePoy & 
Coffee 2001; Suchkov & Schultz 2001). 

5.3 Number of planets detectable via microlensing 

Detailed monitoring of microlensing events has so far failed 
to find anomalies in the light curves characteristic of mas- 
sive, bound planets (Bennett & Rhie 1996; Albrow et al. 
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2001; Gaudi et al. 2002). The results exclude the possibility 
that more than around a third of lensing stars have Jupiter 
mass planets at radii 1.5AU < a < 4AU. Our results are 
consistent with this limit. For masses greater than O.SMj, 
we estimate that 7% of the stars targeted by radial velocity 
surveys have planets within this range of radii. If this is also 
true of the lensing stars (which are typically less massive - 
around 0.3Mq), then it suggests that the existing monitor- 
ing campaigns are within a factor of ~ 5 of reaching limits 
where a detection may be expected. 

5.4 Probability of forming a planet which does 
not migrate 

Unlike in the case of extrasolar planetary systems, the ev- 
idence for migration in the Solar System is extremely lim- 
ited. In situ analysis of the Jovian atmosphere shows that 
it may contain material that originated in the outer Solar 
System (Owen et al. 1999), while substantial migration of 
Uranus and Neptune has been suggested (Thommes, Dun- 
can & Levison 1999). There is little to suggest that Jupiter 
formed a significant distance away from its current location. 

Migration models permit both inward and outward mi- 
gration, so forming a planet which does not migrate is per- 
fectly possible, though it requires fortuitous timing. If we 
take 'no significant migration' to mean that the planet moves 
a distance, 

^<0.1, (11) 

then in our standard disc model the window of formation 
time over which this occurs lasts for At — 6.8 x 10 4 yr. This 
may be compared with the 0.5 Myr range of formation times 
that result in the planet having a semi-major axis in the 
range 0.1AU < a < 6AU. Making the same assumption 
as before - that the rate of planet formation remains con- 
stant across this larger interval - we estimate that no signif- 
icant migration occurs in approximately 10-15% of systems 
in which a planet survives. The no migration outcome which 
may describe Jupiter in the Solar System would then be an 
unlikely event, but not a rare one (Lineweaver & Grether 
2002). Of course forming several planets, none of which mi- 
grated significantly, would be much less likely. However, as 
already noted, migration of the other giant planets in the 
Solar System has been seriously considered and does not 
appear to be excluded on observational grounds. 

5.5 Eccentricity 

Our model for migration ignores the typically substantial ec- 
centricities of extrasolar planets. How serious this omission 
is depends upon the (unknown) mechanism which leads to 
the eccentricity. If interactions between a single planet and 
the disc are responsible (Artymowicz et al. 1991; Papaloizou, 
Nelson & Masset 2001; Goldreich & Sari 2002), then a de- 
scription similar to the one we have developed here ought 
to suffice. Conversely, eccentricity may arise from planet- 
planet scattering in a multiple planet system (Rasio & Ford 
1996; Weidenschilling & Marzari 1996; Lin & Ida 1997; Pa- 
paloizou & Terquem 2001). This would be accompanied by 
order unity changes in the semi-major axis of the surviving 
planet. The predicted distribution of planetary orbital radii, 



shown in Figure g, would be changed if all planets (say) suf- 
fered a shrinkage of their orbits by a factor of two subsequent 
to migration. We have no additional need to invoke such a 
process to explain the orbital radii of the planets, but more 
work in this area is warranted. 



6 CONCLUSIONS 

We have shown that the radial distribution of massive extra- 
solar planets, beyond 0.1 AU, is consistent with the planets 
forming at 5 AU, before migrating inwards through a vis- 
cous protoplanetary disc. This means that although theory 
suggests that massive planets could form closer to the star 
(Papaloizou & Terquem 1999; Bodenheimer, Hubickyj & Lis- 
sauer 2000; Sasselov & Lecar 2000), there is no observational 
requirement for planets to form at a < 3AU. The predicted 
distribution is a robust feature of disc models which include 
dispersal of the gas via a disc wind, so the level of compar- 
ison is currently limited by the small numbers of detected 
planets. The continuation of existing radial velocity surveys, 
along with forthcoming astrometric searches (e.g. Lattanzi 
et al. 2000) , will allow for more stringent tests, provided that 
the selection effects of the surveys are well understood. 

The good agreement between observations and the- 
ory allows us to use the model to estimate the fraction 
of stars that possess massive planets, including those cur- 
rently undetected on account of their low mass or large 
semi-major axis. Limiting ourselves to planets with masses 
0.3Mj < M p sin(i) < 10A/j, and radii 0.1AU < a < 5AU, 
we estimate that around 15% of stars in the current tar- 
get sample of extrasolar planet searches possess planets. In 
some of our models a significant number of planets migrate 
outwards, and these could form a sizeable additional pop- 
ulation. These numbers are comparable to those found by 
Trilling, Lunine & Benz (2002), and imply that even in ex- 
isting samples, most of the massive planets remain to be 
discovered. 

Finally, we have used the results to quantify the fraction 
of systems in which the most massive planet does not suffer 
significant radial migration. This allows us to address the 
concern that migration theories could be inconsistent with 
the detailed knowledge of our own Solar System, in which 
Jupiter may have formed close to its current location. We 
find that a no-migration outcome requires fortuitous timing 
of the epoch of planet formation, but not extraordinary luck. 
Of the order of 10-15% of planetary systems are expected to 
have a massive planet which has migrated by less than 10% 
in orbital radius. Having a most massive planet at the orbital 
radius of Jupiter is thus expected to occur in around 1-2% 
of solar-type stars. The orbital radius of Jupiter, at least, 
is consistent with our understanding of extrasolar planets, 
though the near-circular orbit may yet prove to be unusual. 
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